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ABSTRACT 

The equation of state of hydrogen, helium, and water effects interior structure 
models of giant planets significantly. We present a new equation of state data 
table, LM-REOS, generated by large scale quantum molecular dynamics simula- 
tions for hydrogen, helium, and water in the warm dense matter regime, i.e. for 
megabar pressures and temperatures of several thousand Kelvin, and by advanced 
chemical methods in the complementary regions. The influence of LM-REOS on 
the structure of Jupiter is investigated and compared with state-of-the-art results 
within a standard three-layer model consistent with astrophysical observations of 
Jupiter. Our new Jupiter models predict an important impact of mixing effects of 
helium in hydrogen with respect to an altered compressibility and immiscibility. 

Subject headings: planets and satellites: individual (Jupiter) - equation of state 

1. Introduction 

Jupiter consists by more than 85% in mass of hydrogen and helium. Among the minor 
constituents, oxygen contributes the largest fraction. Although the major constituents are 
rather simple elements, they occur almost exclusively in extreme, up to now largely unknown 
states of matter. For instance, about 90% of the planetary material in Jupiter is in a 
high-pressure state beyond the 1 Mbar level up to about 70 Mbar in the center, and the 
temperature varies by two orders of magnitude from the cold and dilute outer envelope up 
to about 20.000 K in the deep interior. In terms of usual plasma parameters, the interior 
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covers the transition from weak to strong ion coupling, and from weak to moderate electron 
degeneracy without relativistic effects. Giant planets like Jupiter are therefore particularly 
suitable to study warm dense matter (WDM), especially mixtures of H, He, and H2O. 

In this paper we apply for the first time ab initio equation of state (EOS) data for 
H, He, and H2O, representing at least 97% of the planetary material in Jupiter, to interior 
models in order to derive more reliable implications for an improved structure model as well 
as for the H/He phase diagram. We have performed large scale quantum molecular dynamics 
simulations for the WDM region and applied advanced chemical models in the complemen- 
tary regions. These EOS data for the major components were combined by linear mixing 
into a new data table LM-REOS (Linear Mixing Rostock Equation of State). Fundamental 
problems concerning hydrogen and helium such as the location of the nonmetal-to-metal 
transition (NMT), whether or not this transition is accompanied by a thermodynamic phase 
instability - the plasma phase transition (PPT) -, and the existence of a miscibility gap are 
at the same time key issues for the construction of planetary models. 

Strongly improved high-pressure data for planetary materials as well as new and more 
accurate observational data for Jupiter and the other planets have motivated intensive efforts 
to model giant planets. We outline briefly major steps in the development of planetary 
models that lead to our present understanding of Jupiter's interior structure. 



In one of the first approaches, iHubbardI (119681 ) started from estimates of Jupiter's net 



heat flux and inferred a slightly superadiabatic temperature gradient, causing convective 
instability throughout the planet. With his fully homogeneous and, thus, overdetermined 
model he compared with cu rrent observatio nal constraints in order to discriminate between 
different equations of state ( Hubbardlll969 ). 



Stevenson fc Salpeter ( 1977a| jbl) examined the transport properties of hydrogen- helium 



mixtures and the influence of a PPT and of H/He phase separation on the structure and 
evolution of hydrogen-helium planets. Either of these processes would divide the fluid interior 
into two convective layers, a molecular outer envolope and a metallic inner envelope with 
different helium abundances. The NMT was predicted to occur around 2-4 Mbar, and helium 
immiscibility between 4-40 Mbar and below 8000 K. Depending on the true position of the 
NMT coexistence line and the demixing curve in the H/He phase diagram, the extension 
of the inhomogeneous, diffusive transition region in between could be very narrow. Minor 
constituents were argued to be redistributed in consequence of helium immiscibihty. 

Voyager data indeed unveiled a depletion of helium in the outer molecular envelope, 
indicating that H/He pha se separation has happened there. A thermodynamically consistent 



H-EOS was developed by lSaumon fc Chabrierl (119921 ) for astrophysical applications within a 
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chemi cal model that predicted a discont i nuous NMT. This PPT has been proposed for a long 

igram 



time ( 


Norman & Starostin 


1968 




19' 




Zeldovich & Landau 


1944 


) and the phase diag 


was discussed in detail ( 


Beule et al. 


2001 




Ebeline & Norman 


2003: 


Hoist et al. 


2007) 



first experimental signature for a PPT in deuterium was given quite recently (IFortov et al. 



20071 ). Using this H-EOS, the observational constraints for Jupiter could only be adjusted 
if the envelopes were allowed to have a de nsity discontinuity, a nd if a several Earth masses 
sized core of rocks and ices was admitted ( Chabrier et al.l [19921 ). For simplicity, helium was 
(and in recent models still is) added to the EOS of hydrogen via the additive volume rule, 
and the fraction of helium in both envelopes was used to adjust the gravitational moments 
J2 and J4. Other constituents than hydrogen and helium - summarized as heavy elements 
or metals - were represented by helium too. 



A significant step forward in constructing r eliable Jupiter naodels was done by iGuillot 
f|l999al lb[) who used an improved co mput er code (IGuillot fc Morellll995l ). improved EOS data 
of Saumon, Chabrier & Van Horn (119951 ) (SCvH), and more ac curate observational d ata for 
the atmospheric helium abundance from the Galileo mission (jvon Zahn et al.lll998l ). This 
improved Jupiter model predicted a core mass of 0-10 Earth masses (Mq), a total amount 
of heavy elements of 11-41 M®, and a heavy element enrichment of up to 6.5 times the solar 
value in the molecular region. The great margins resulted from a variation of the observables 
within their error bars, but most of all from the uncertainty in the H-EOS data itself in 
the region around the NMT. By choosing appropriate pair potentials for the H2 molecules, 
the underlying H-EOS (^H-SCvH-p pt) was constructed to reproduce sho ck compression data 
derived from gas gun experiments (INellis et al.l Il983l : iThiel et al.l Il973l ) . A first-order phase 
transition was found around 2 Mbar, the slightly decreasing coexistence line ends in a second 
critical point at about 15000 K. Its interpolated version (H-SCyH-i) has been widely use d 
in planeta ry modelling ([Gudkova fc Zharkov 1999 : Guillot 1999a : Saumon fc Guillot 2004). 
As before ( IChabrier et al.lll992[ ). the fraction of metals was obtained from the excess density 
with respect to that of a given H/He mixture. 



In an alternative approach, iKerleyl (l2004al ) calculated three-l ayer models ( a core of rocks 
or ices an d two fluid env elopes) using a revised Sesame EOS for H (lKerleyll2003l ) . an improved 
He-EOS jKerlevll2004b[ ). and the PANDA code for the icy components H2O, NH3, CH4, and 
H2S and their atomic constituents to represent r netals. Choosing threefold enrich ment in 
heavy elements in accordance with measurements (lAtreya et al.ll2003l : lMahaffyll2000t ) for the 
molecular envelope, an enrichment of 7.5 times solar was required in the metallic envelope 
to match J2, in total 35 M(^. The large heavy element enrichment in the inner envelope was 
argued to result from a dispersion of planetesimals during the formation process. 



The latest extensive calculations of Jupiter models have been performed (jSaumon &: Guillot 
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20041 ) using the Sesame EOS 7154 of water and the Sesame EOS 7100 of dry sand to represent 
metal s within the envelop es and for hydrogen using the simple linear mixing (LM) model of 
Ross (IHolmes et al.lll995l ). The LM model served as a basis for the generation of different 
H-EOS representing the current uncertainty on the experimental deuterium EOS data. Not 
all of the H-EOS applied gave satisfactory results. While acceptable Jupiter models were 
placed within the range of the previous models cited above, some LM-variants as well as 
the original Sesame EOS 5263 were not compatible with Jupiter's gravitational properties 
or age. 



Other attempts (IHolst et al.ll2007l ) to construct Jupiter models with the Ses ame EOS of 
hydro gen or a combination of fluid variational theory with plasma effects (EVT"*") (IHolst et al 
20071 ) also failed to reproduce the observational constraints. 



To decide whether a mismatch between calculated and observed parameters is due to an 
inaccuracy of the EOS or due to improper assumptions in the planetary model, accurate ex- 
perimental data for the high-pressure phase diagram of hydrogen and hydrogen-helium mix- 
tures are much-needed. Various experiments have been performed in the recen t years to study 



the complex physical behavior of hydrogen in the region up to several Mbar (IBoriskov et al. 



20051 : ICoUins et al.lll998l : iKnudson et al.ll2004j : IWeir et al.lll996l ). The transition from a non- 
conducting, r aolecular fluid to an atomic, conducting flui d can be deri ved from electrical 
conductivity (INellis et al.lll996l ) as well as reflectivity data (ICelliersll2000l ). For instance, the 
NMT in fluid hydrogen at temperat ures of about 300 K has been pinpointed at 1.4 Mbar 
in multiple shock wave experiments (jWeir et al.lll996l ). Reflectivity rn easurements indicate 
an NMT along the Hugoniot curve at about 0.5 Mbar (jCelliersll2000l ). Most recently, the 
flrst experimental signature of a flrst-order phase transition connected with t his NMT has 



been found in quasi-isentropic shock compression experiments in deuterium (IFortov et al. 



20071 ). The maximum compression ratio along the principal Hugoniot curve of hydrogen 
is not u ltimately flxed but a value of 4.25 has be en found in various shock- wave exper- 
iments teoriskov et al. 2005; Knudson et al. 2004) which is not predicted by any of the 
theoretical H-EOS mentioned above except the Sesame EOS. 

On the other hand, current ab initio calculations yield not only agreement with experi- 



200d:lDesiarl3|2003yBonev et 



2001 



ment al Hugoniot curves (iLenosky et al 
20081) but also with reflectivity data (jCoUins et al. 

20081 ) at high pressures. Therefore, we use in this study ab initio EOS data derived from 



Mazevet et al 



2004 



200 



i 



Hoist et al. 



Hoist et al. 



extensive quantum molecular dynami cs (QMD ) simulations for the mo st abundant planetary 
mate rials hydrogen ( Hoist et al.ll2008 ). helium ( Kietzmann et alll2007 ). and water ( French et al, 
20081 ) and study the impact of these data on models of Jupiter within the three-layer struc- 
ture assumption. 
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Our paper is organized as follows. In the following section we describe our procedure to 
calculate three-layer planetary models. In §3] the new LM-REOS data tables for H, He, and 
H2O are introduced and the differences to the widely used SCvH EOS are exemplified by 
means of the Jupiter isentrope and the relevant plasma parameters. Our results for Jupiter 
are presented in ^ First, our numerical results are checked by using the SCvH EOS. We then 
process the new LM- REOS data tab l es and compare th e predictions for Jupiter's interior 
with previous results (iGuillotI Il999al : iGuillot et al.l 120031 ) and with additional calculations 
using the SCvH-ppt EOS regarding the distribution of chemical species. We discuss our 
results for Jupiter in ^ and propose key issues for future work on high-pressure EOS data 
and improvements of planetary models. A summary is given in ^ 



2. Three-layer planetary models 
2.1. Parameters for Jupiter 



Our interior models follow the approach described in detail by iGuillotI (jl999al ). The 
planet is assumed to consist of three layers: two fluid, homogeneous, and isentropic envelopes 
composed of hydrogen, helium and metals and an isothermal core of rocks or ices. 

For given EOS of these materials, the precise fractions of helium and metals within the 
envelopes and the size of the core are constrained by the requirement of matching the ob- 
servational constraints, see Tab. [TJ Jupiter's mass-radius relation (total mass M.t, equatoria l 



radius Req) requires the dominant mass fraction to be attributed to hydrogen (lGuillotll2005l ). 
The mean helium fraction Y is assumed to reflect the hydrogen-helium ratio of the proto- 
solar cloud. Due to the assumption of homogeneous layers the atmospheric helium fraction 
Ymoi is assigned to the whole outer envelope, starting at the 1 bar level with temperature 
Ti. The adjustment of the gravitational moments J2 and J4 assuming rigid rotation with 
angular velocity u sets limits to the fraction of metals Z^^i in the outer envelope and Z^^t 
in the inner envelope. Observational error bars translate into uncertainties of the resulting 
values of these structure properties, allowing for a variety of acceptable models. Compar- 
ing acceptable models with other predictions, such as from formation theory or abundance 
measurements, allows to evaluate the models, see §3 
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Table 1. Measured data that are considered in the modelhng procedure. 



Parameter 



Jupiter 



mass Mj [g] 
equatorial radius Reg [m] 
average helium mass fraction Y 
helium mass fraction Y^^i 
temperature Ti [K] 
angular velocity u 
J2/IO-2 



7 b 



1.8986112(15) X 10 
7.1492(4) X 10 

n 97c:+0.003 c 
U.z ( O_o.o05 

0.238(7) <^ 
165-170 
27r/9h55 ^ 
1.4697(1) ^ 
-5.84(5) ^ 



30 



^Numbers in parenthesis are the uncertainty in the last digits of the given value. 



^taken from lOuillot et all (120031 ). 



teahcalll (119951). 



von Zahn et all (119981 ). 



"taken from [GuiUotl (Il999ah 
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Besides being the simplest interior structure model that is consistent with the observa- 
tional constraints mentioned above, the three-layer stucture has been discu ssed and evalu - 
ated in connection with various underlying or neglected physical properties (lGuillotlll999al). 
For in stance, a small opacity could occur in a narrow region around 1500 K (iGuillot et al. 
1994al jb[). causing a stable radiative layer and thus lengthen the cooling time of the planet 
without affecting the element distribution. Differential rotat ion is observed in the atmo- 
sphere and could modify J2 up to 0.5% and J4 up to 1% (IZharkov fc Trubitsynl Il978l ). 
Assuming more than two envelopes would probably reflect the real structure of Jupiter more 
adequately, but also entail more degrees of freedom than can be adjusted by the constraints 
presently available. 

Layer boundaries that devide the fluid planet into convective layers with different com- 
positions have been predicted theoretically from general thermodynamic prope rties of warm 
dense hydrogen and hydrogen-helium mixtures ( Stevenson fc Salpeter 1977a| jbl). Possible 
sources of inhomogeneity from thermodynamic considerations are phase separation due to 
first-order phase transitions or immiscibility between different components, and from evo- 
lution theory the presence of a slowly eroding core. In case of Jupiter, the PPT in hy- 
drogen is a candidate for an important first-order phase transition, besides those for the 
minor constituents as, e.g., H2O in the outermost shells. Phase separa tion may occur be- 



tween neutral helium an d metallic hydrogen near the PPT of hydrogen (IKlepeis et al.lll991 



Pfaffenzeller et al.lll995l ). Both effects motivate a three-layer structure assumption. 



2.2. Basic equations 



The construction of acceptable interior models is performed in an iterative procedure. 
Each model is defined by a set of distinct parameter values to be matched. For the accurately 
known Mj, Req,uj, J2, these values are the observed ones within an error of one standard 
deviation (la); for the less accurately known Y and J4, the values are chosen in the la 
range with a computational error of less than 10% of their la error. The interior models 
are not optimized to meet Jq. Resulting values of Jq lie always in the la range, if J2 and J4 
do. For a H-EOS with PPT, the transition pressure Pm from the molecular to the metallic 
layer is chosen to coincide with the entry of the outer isentrope into the instability region. 
For a H-EOS without PPT, the pressure is varied between 1.4 and 5 Mbar. For a distinct 
choice of mass fractions X of hydrogen, Y of helium, Z of metals, the EOS of the mixture 
is calculated, see §3.61 From these tabulated thermal and caloric EOS of the mixture, the 
entropy can be calculated via thermodynamic relations. The isentropes of the outer and 
of the inner envelope and the isotherm of the core are defined by the boundary conditions 
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(Ti,Pi), (Tni,Pm), and (Tc,Pc), where the transition temperature Tm results from the outer 
isentrope at P = Pm, and the core temperature from the inner isentrope at P = Pc. 
Since the size of the core is not known in advance, Pc is varied within the iteration scheme 
described below until the total mass condition holds. Along this curve of piecewise constant 
entropy or temperature the equation of motion describing hydrostatic equilibrium, 

p-ivp = v(v + g), (1) 

is integrated inwards, where p is the mass density and U = V + Q the total potential 
composed of the gravitational potential 

nO=G/<iVi^ (2) 

and the centrifugal potential Q defined by 

- VQ(f) =ujx {Oxr) . (3) 

A proper description for the potential of axisymmetric but, due to rotation, oblate planets is 
a multipole expansion, which transforms the general expressions ([2]) and ([3]) into the following 
forms with spherical coordinates r,6, and Legendre polynomials P2n{t = cos 6*) 

V{r, e) = -J2 P2n{t) / d'r'pir', 6') (-] PUt') , (4) 

^ n=0 V ^ / 

Q{r,e)=^-uh\l-P2{t)) . (5) 

The index k = n describes the external gravitational potential for r > r' and k = —{n + 1) 
the internal gravitational potential for r < r'. At the surface, the internal gravitational 
potential vanishes and V reduces to 

V{r, ^) = ^ (l - E (^) " J2nP2n{t)^ , (6) 



with the gravitational moments 



"eg 



J2n = j^ ld'r'p{r\ey'-P,^{t'). (7) 



Calculating the integrals in equations (jl]) and ([7]) requires know 
and the density distribution. We apply t he Theory of Figures ( 



edge of the planetary shape 



Zharkov fc Trubitsynlll978l ) 



and the method described therein and in (INettelmann et al.l 120071 ) for solving the system of 
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integral equations on computers to iteratively calculate the planetary shape, the potential, 
and the density distribution. The key point of this theory is to reduce the dimension of the 
problem by replacing the former radial coordinate r by its representation 

(oo 
l + 5^S2.(/)P2.(t) 
i=0 

on equipotential surfaces (/ =const) with expansion coefficients S2i{l)- Reasonable and prac- 
ti cal cutof f indice s for the sum in equation ([8]) are 3, 4 or 5. In accordance with the treatment 
of lGuillotI (Il999al ). we consider only coefficients up to third order, i.e. i = 3, which introduces 
an error into the calculation of the gravitational moments of the order of magnitude of Je ot 
about five times the observational error of J2. In terms of the level coordinate I = . . . R, 
the equations to be integrated along the piecewise isentrope take the simple form 

1 dP{l) dU{l) 



p{l) dl dl 

for the equation of motion and 



(9) 



^ = 4vr/V(0 (10) 

for the equation of mass conservation. The mean radius R of the surface (at the 1-bar 
level) is not known in advance and has to be adjusted in order to satisfy r{6 = tt/2) = 
Reg. The coefficients S2i in equ ation (IHl), the figure function s, determine the shape of the 
rotating planet. As suggested by lZharkov fc Trubitsynl (119781 ) they are calculated for a fixed 
density distribution p{l) in an iterative procedure until convergence of at least 0.01% for 
S2 is reached, which turns out to be accompanied by an error of 0.05% for S4 and 0.1% 
for Sq. The converged figure functions enter into equations (Q, (fTOj) . and ([7]) and a revised 
density distribution is calculated resulting in a different core mass and, thus, a different 
central pressure Pc. After about six iterations of calculating converged figure functions and 
corresponding density distributions, J2 and J4 are converged within their observational error 
bars, but not necessarily to the observed values ^2^*^, Jf^^- To achieve agreement also with 
Y, the guess of Zmoi, ^mct and Fmet is improved slightly and the procedure decribed above is 
repeated, starting with the present density distribution, until the desired values of J2, J4 and 
Y have been approached. This procedure can, in principle, be applied to all giant planets in 
the Solar system. 



2.3. Accuracy 



In order to give reliable results for the interior structure model and to study the infiuence 
of different EOS data sets, the numerical treatment has to ensure a definite precision. The 
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accuracy of acceptable interior models should be in the order of or smaller than the relative 
error of the most accurate observed quantity to be adjusted during the procedure. For Mj, 
Reg and J2 the observational error is below 0.01% and for F, Fmob ^"1, and J4 larger than 
1%. The uncertainty associated with the QMD data themselves induces an error of the order 
of 1% to the isentropes. Hence we consider a numerical accuracy of 0.1% as sufficient to 
study the effect of the uncertainty of the less accurately known observables and for conclusive 
interior models regarding the EOS applied. 

Numerical errors occur due to the integration of equations (I9l) and ( ITOl) . the differen- 
tiation of the gravitational potential, and the integrals in equation (jlj), transformed to the 
level coordinate /. We have tested the accuracy of integrating the differential equations ([9]) 
and flTUl) for a non- rotating polytrope of index 1. Our numerical results for the profiles m(r) 
and P(r) obtained via a fourth-order Runge-Kutta method with adaptive stepsize control 
in regions with steep gradients differ by less than 0.001% from the analytical solution in the 
outer 80% of the planet. In the inner 20% of the planet including the core region, a shooting 
splitting method is applied to localize the outer boundary of the core. We then search for a 
solution starting from the center to meet the envelope solution there with a deviation of less 
than 0.1%. Since a polytropic model does not have a core, this method cannot be applied 
there and the difference rises to 50% for the mass and 0.1% for the pressure near the center. 

To estimate the error resulting from the calculation of the integrals and its enhancement 
during the iteration procedure, we have compared fully converged planetary models. The 
resulting mass fractions Z^o\ and Zmet and the size of the core differ by less than 2% if 
the number of intervals dV in the integrations is doubled. Since density discontinuities at 
layer boundaries are extended over a small but finite intervall, they do not cause numerical 
difficulties. Extending this intervall from 0.001 i?j to 0.01 i?j effects the resulting core mass 
by about 0.1 and the heavy element abundance by about 2%. 

We conclude that a convergence of our numerical procedure within 0.1% can be ensured, 
but the resulting values for the model parameters are uncertain within ±0.1 or 2%. 

3. EOS data 

Interior models of Jupiter's present state require EOS data for H, He, and metals from 
about 160 K at 1 bar up to 25000 K at about 45 Mbar, and for core materials (rocks, ices) 
of 35-80 Mbar around 20000 K. No state-of-the-art thermodynamic model is capable of cal- 
culating accurate EOS data for all planetary materials in this wide pressure-temperature 
range. However, various EOS data sets have been constructed for this purpose by adopting 
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known limiting cases and interpolating in between or extending them till phase bound- 
aries. Examples are the Se same EOS data tables 5251, 5761 and 7150 for H, He and 



H2O JSESAME hbrarvl 



19921) and c hemical models such as the SCvH EOS (ISaumon et al. 



I995I ) and FVT+ (IHolst et al.l 120071 ). These data sets have the largest uncertainties in the 
WDM regime. We have, therefore, performed extensive QMD simulations to generate accu- 
rate EOS data for warm dense H, He and H2O, i.e. for moderate temperatures and densities 
higher than 0.2 (H, He) or 1 g cm~^ (H2O). The new ab initio data cover at least 97% of 
the planetary material inside Jupiter. For lower densities, the simulation times increases 
dramatically so that other methods have to be applied. 

We outline the method of QMD simulations below and describe the construction of 
improved EOS tables for H, He, and H2O by combining QMD results for WDM with chemical 
models in complementary regions of the phase diagram. Predictions of these new EOS tables 
are compared with SCvH EOS with respect to the Jupiter adiabat and plasma parameters. 



3.1. QMD simulations 



For the QMD sin aulations we have used th e cod e VASP (Vienna Ab Initio Si mulation 
Package) developed by Kresse fc Hafner f 1993al jbl) and Kresse fc Furthmiiller ( 1996 ). Within 
this method, the ions are treated by classical molecular dynamics simulations. For a given 
ion distribution, the density distribution of the electrons is determined using finite temper- 
ature density functional theory (FT-DFT). The electronic wave functions are represented 
by plane waves and calculated using projector augraented wave (PAW) pseudopotentials 



(]Blochllll994J : iKresse k Joubertlll999l : iDesiarlais et all 120021 1. The central input into DFT is 
the exchange-correlation functional accounting for interactions between electrons and ions 
and the electrons themselves. It is calculated wit hin g eneralized gradient approximation 



(GGA) using the parametrization of iPerdew et al.l (119991 ). 



The convergence of the thermodynamic quantities in QMD simulati ons is an important 
i ssue. We found that a plane wave cut off energy o f 700 ( He), 900 (H2O, see lMattsson fc Desjarlais 
( l2006l )). and up to 1200 eV (H, see iDesjarlaid (120031 )) is necessary to converge the pres- 
sure within 3% accuracy. Furthermore, we have checked the convergence with respect to 
a systematic enlargement of the k-point set in the representation of the Brillouin zone. 
Higher-order k points modify the EOS data only within 1% relative to a one-point result. 
Therefore, we have r estricted our calcu lations to the F point for water or the mean value 
point (1/4, 1/4, 1/4) (IBaldereschil Il973l ) for H and He. The overall uncertainty of the QMD 
data due to statistical and systematic errors is below 5%. The simulations were performed 
for a canonical ensemble where the temperature, the volume of the simulation box, and the 
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particle number in the box are conserved quantities. We consider 32 to 162 atoms in the 
simulation box and perio dic boundar y conditions. The ion temperature Tj is regulated by a 
Nose- Hoover thermostat (iNosg Il984l ) and the electron temperature (wi th Tp = Tj) is fixed 
by Fermi weighting the occupation of bands using Mermins approach (IMerminl Il965l ). At 
high temperatures, the number of bands to be considered increases exponentially with the 
electron temperature so that a treatment of high-temperature plasmas is rather challenging. 
However, simulations of WDM states at several 1000 K as typical for planetary interiors are 
performed without any serious difficulties. When the simulation has reached thermodynamic 
equilibrium, the subsequent 500 to 2000 time steps are taken to calculate the EOS data as 
running averages. 



3.2. Construction of the hydrogen EOS 



For densities below 0.1 g cm~^ we have used FVT+ which treats the dissociation of 
molecules self-consistently and ta kes into account the ionization of the hydrogen atoms 
via a respective mass actio n law (iHolst et al.l 120071 ) . F VT"*" makes use of Fade formulas 
(IChabrier &: Potekhinl Il998l ) for the plasma properties. Below 3000 K ionization does not 
play a role at planetary pressures and is thus neglected. FVT+ data on the one hand and 
QMD data on the other hand are then combined within a spline interpolation. 

Both methods yield different reference values Uq for the internal energy, which have 
been fixed by the Hugoniot starting point of fiuid, molecular hydrogen at 19.6 K and 0.0855 
g cm~^. The energy difference is added to the QMD data in order to get a smooth transition 
between the data sets. The final hydrogen EOS is referred to as H-REOS. Examples of 
isotherms are shown in Fig. [H 



3.3. Construction of the helium EOS 



The He-EOS named He-REOS is a combination of QMD data for densities be tween 0.16 

and 10 g cm~^ and temperatures between 4000 and 31600 K, and Sesame 5761 data (ISESAME library 
19921 ) for densities below 0.01 and above 16 g cm~^ and for lower and higher temperatures. 
For intermediate densities, the isotherms of both data sets were interpolated. In Fig. [2] we 
compare isotherms of He-REOS with Sesame 5761. QMD data for helium indicate a slightly 
lower compressibility around 1 Mbar. Optical properties, especially the metallization of fiuid 
helium, as relevant for giant planets and white dwar f atmospheres have been investigated 
recently (IKietzmann et al.l 120071 : iKowalski et al.ll2007l ) by means of QMD simulations. 
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3.4. Construction of the water EOS 



The EOS of water named H2O-REOS is a combination of four EOS data sets. QMD 
data are considered for 1000 <T < 10000 K and between 2 and 7 g cm~'^, as well as for 10000 
< T < 24000 K and between 5 and 15 g cm~^. Details of the QMD simula tions for H2O and 



the respective EOS data will be published elsewhere ( French et al. 2008[). For the phases 



ice I and liquid water we apply accurate data tables named FW (IFeistel fc Wagned l2006l ) 
and WP (IWagner fc Prusa 120021). respectively. All other regions are filled up with Sesame 
7150 data (jSESAME librarylll992l ). Intermediate EOS data are gained via interpolation. 



Fig. El shows six pressure isotherms for H2O-REOS representing different phases of H2O 
and the pressure-density relation of the water fraction representing heavy elements along a 
typical Jupiter isentrope. Due to the phase diagram of pure H2O, in the outermost region 
of Jupiter the H20-component transits from ice I to liquid water below 300 K and further 
to the vapor phase around 550 K. At about 4000 K, H2O again reaches densities of 1 gcm~^ 
along the Jupiter isentrope. Above this poi nt, water shows strong mole cular dissociation 
and becomes also electronically conductive (IMattsson fc Desjarlaid l2006l ). From here on, 
the EOS of the water component in Jupiter is described by our QMD data. Following the 
isentrope further into the deep interior, water remains in the fluid, plasma-like phase. 



Our H2O phase diagram smoothly adds to recent QMD simulations by lMattsson fc Desjarlais 

( I2OO6I ). The EOS data differ at densities between 3-10 g cm~^ from the pure Sesame EOS. 



3.5. EOS of metals and core material 



We make use of the EOS table for water to represent metals in the envelopes. Alterna- 
tively, the EOS table of helium scaled in density by a factor of 4 (He4) is used to represent 
the atomic weight of oxygen. Oxy gen is the third abundant element in the sun and probably 
also in Jupiter, see iKerleyl (l2004al ) for a discussion. We compare the respective isotherms of 
the EOS for metals in Fig. H] in that density and temperature region where QMD data of 
H2O are applied. The H2O pressures are up to 10 times higher than the He4 data, probably 
resulting from enhanced ideal and pronounced Coulomb contributions in water plasma. 

Previous models of J upiter assumed that the core consists of rocks surrounded by 
ices ( IChabrier et al.lll992l : iGudkova fc Zharkovl Il999l : IStevensonI Il982l ). Traditionally, rock 
material is a mixture of S, Si, Mg, and Fe while ice refers to a mixture of C, N, O, and H. 



Due t o the core accretion scenario for the formation process of giant planets (lAlibert et al. 



2OO4J : iPoUack et al.lll996l ). these elements have mainly entered the planetary interior in form 
of icy compounds of CH4, NH3, H2O, and of porous material composed of Si02, MgO, FeS, 
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and FeO. Their present state in the deep interior of the planet cannot be inferred from ob- 
servations or interior model calculations. Our QMD simulations show clearly t hat H^O is in 



a pla sma and not in a solid phase at conditions typical for the core of Jupiter (jPrench et al. 



20081 ). Since CH4 and NH3 are more volatile than H2O, they are expected to prefer the 



plasma phase too in the core region of Jupiter. Because a smaller-sized core results from 
interior models r ather than frorn formation models of Ju piter, the core has been argued to 
erode with time (lAlibert et al.ll2004 iGuillot et al.ll2003l ). Motivated by the possibility of 
core erosion and by the miscibility of the icy component within the metallic envelope we 
assume that the core consists purely of rocks in all our calculations for Jupiter. The EOS of 
rocks is adapted from lHubbard k. Marley Jl989h . 



3.6. Mixtures 



The EOS of a mixture with mass fractions X of hydrogen, Y of helium, and Z of metals 
i s calculated by means of the EOS for the pure components using the additive volume rule 
(jChabrier et al.l Il992l : [Peebles! Il964j ) for the internal energy. 



n(P, T) = Xux{P, T) + Yuy{P, T) + Zuz{P, T) 



and the mass density. 



X 



+ 



Y 



+ 



p(P,T) px{P.T) Py{P.T) pziP,T) 



(12) 



on a grid of pressures P and temperatures T. Non-ideality effects of mixing are thus not 
taken into account. 

In this paper we refer to three different linear mixtures of EOS data abbreviated by LM- 
REOS, SCvH-ppt and SCvH-i. The new data set LM-REOS consists of a linear mixture as 
described above of H-REOS, He-REOS, and of H2O-REOS or He4-RE0S for metals. SCvH- 
ppt as used in this work consists of H-SCvH-ppt for hydrogen and He-SCvH for helium 
( ISaumon et al.lll995l ). and of He-SCvH scaled in density by a factor four as representative 



for m etals. SCvH-i as used in this work consists of the interpolated version (jSaumon et al. 



I995I ) of H-SCvH-ppt together with the same EOS for He and metals as for SCvH-ppt. Apart 
from SCvH-ppt, none of the hydrogen or helium EOS used here shows a first-order phase 
transition. 

For given mixtures, the piecewise isentrope defined by the boundary condition is cal- 
culated as described in §2.21 When using H2O for metals we make the following exception 
for temperatures below 1000 K. At these small temperatures where H2O passes the phases 
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ice I, liquid water, and vapor, we first calculate the isentrope of the H/He mixture and then 
add the desired mass fraction of H2O to the -P — p relation of that isentrope according to its 
P — T relation. At 1000 K, the pieces of isentropes for smaller and higher temperatures add 
very smoothly, as can be seen from Fig. [31 



3.7. Comparison with chemical models 

Results of our LM-REOS for the coupling parameter F, the degeneracy parameter G, 
and the pressure-density relation along the hydrogen adiabat are compared with the chemical 
models SCvH-ppt and SCvH-i. 



We have applied the definition of F and G given by IChabrier fc Potekhiru (119981 ) for 



multi-component plasmas with mean ion charge (Z) and electron density Ue 



Ai^eokB \ 3 

G = ||2m,(3.)V3,^_ (14) 

These plasma parameters are calculated along typical Jupiter adiabats by taking into ac- 
count the ionization degree of hydrogen and helium. For models with the SCvH-ppt and 
SCvH-i EOS, the fraction of H+, He^, and He^^ was taken simply from their tables. This 
is not directly possible for models using the QMD data sets because the underlying strict 
physical picture does not discriminate between bound and free states. However, we have 
estimate d the ionization degree of the He component by using the COMPTRA04 program 



package (iKuhlbrodt et al.ll2005l ). In the underlying chemical model of COMPTRA04, the 
composition of a partially ionized plasma at a given temperature and mass density is de- 
termined by the solution of coupled mass action laws including non-ideal contributions to 
the chemical potential via Fade-formulas and a polarization potential accounting for the 
interaction between electrons and neutral atoms. For the hydrogen subsystem, the degree 
of dissociation is evaluated via the proton-proton distribution function and the coordination 



number derived from the QMD runs; this method is described in more detail by iHolst et al. 



( I2OO8I ). According to state-of-the-art chemical models for hydrogen such as SCvH or FVT+, 
the fraction of atoms is always small so that we can estimate the ionization degree in the 
hydrogen component by taking the dissociation degree as an upper limit. Results for F and 
G are shown in Fig. as function of the mass coordinate in Jupiter. 

At least 80% of Jupiter's mass is in a state with strong coupling and degeneracy effects, 
defined by F > 1 and G < 1. With SCvH-ppt, the step-like increase of F coincides with the 
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density discontinuity at the layer boundary. LM-REOS has a much earher onset of ionization 
in the hydrogen component than both SCvH-ppt and SCvH-i. In the very deep interior, only 
about 5% of the He atoms are singly ionized and less than 1% doubly ionized as inferred from 
the chemical model COMPTRA04 applied to our LM-REOS data, while complete ionization 
is predicted by the SCvH-He EOS so that T becomes shghtly larger in that models. This 
large discrepancy between chemical models in the region of pressure ionization is due to a 
different treatment of correlation effects occuring in the exponents of the mass-action laws 
for ionization, from which the composition is derived. Higher temperatures as with the 
SCvH-ppt EOS diminish the degeneracy and in spite of fully ionized helium the same values 
are found as with LM-REOS. 

More than 50% of Jupiter's mass is in a high-pressure state between 0.1 and 10 Mbar 
where the differences in the hydrogen isentropes are most pronounced and manifest them- 
selves in the structure models, see §4.21 With the transition from FVT^ data to QMD data 
around 0.05 Mbar, the Jupiter isentrope of pure hydrogen becomes much softer than H- 
SCvH-i, with a maximum difference of 30% at 1 Mbar, as shown in Fig. [61 In the innermost 
20% of Jupiter's mass, i.e. above 20 Mbar, the H-REOS isentrope follows the H-SCvH-ppt 
one with 5% higher pressures than SCvH-i. Differences below 0.1 Mbar do not significantly 
contribute to the structure models since this region occupies less than 2% of Jupiter's total 
mass. 



4. Results for Jupiter 

In this section we present results for Jupiter's core mass, for the abundances of metals 
in the envelopes, and for the profiles of the main components along the radius. The size of 
the core Mc determined by the condition of total mass conservation, and the abundances of 
metals Z^oi, ^met which are needed to reproduce J2 and J4 are very sensitive with respect 
to the choice of the EOS and the precise values of t he observation a l para meters to be repro- 



duced. First, we compare our results with those of iGuillot et al.l (120031 ) by using the same 
EOS data sets (SCvH-ppt and SCvH-i) in order to demonstrate that our procedure is sound. 
In §4.21 we then present the range of acceptable Jupiter models using our new LM-REOS 
in comparison with SCvH-based models, explore in §4.31 the behavior of LM-REOS based 
models by varying Ti, J4, and Pm, and discuss in §4.41 the fractions of chemical species calcu- 
lated with LM-REOS and SCvH-ppt for Jupiter models which satisfy the same observational 
constraints given in Tab. [TJ 
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4.1. Consistency with former results 



We follow the notation of iGuillot et al.l (120031 ) who has performed extensive calculations 
for planetary interiors and study four of the models he has introduced explicitly. These 
reference models (A, B, D, E)-TG are based on a three-layer structure as also adapted here. 
Tab. [2] contains their first independent confirmation. In our recalculated models labeled 
(A,. . . , E)-this, we do not allow for a deviation of Y and J4 of more than 1/lOcr; but in our 
models labeled (A',. . . , E')-this, Y and J4 are varied within la in order to achieve the best 
agreement as possible. Corresponding numbers are listed in Tab. [31 
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Table 2. Comparison of present (this) and former (TG) results. 



Model 


H-EOS 


Ti 


Me 


Mz 


Zraet 






SCvH- 


[K] 


[Me] 


[Me] 


[Ze] 


[^o] 


A-TG 




165 


4.2 


33.1 


4.7 


5.3 


A-this 




165 


0.7 


35.6 


6.0 


4.5 


A'-this 




165 


2.8 


34.2 


5.2 


5.0 


B-TG 




170 





35.3 


6.1 


4.1 


B-this 




170 


0.8 


36.4 


6.0 


4.9 


B'-this 




170 





35.5 


6.0 


4.8 


E-TG 


ppt 


165 


4.3 


19.4 


2.5 


2.3 


E-this 


ppt 


165 


6.7 


24.8 


2.9 


3.9 


E'-this 


ppt 


165 


4.6 


22.7 


2.9 


3.3 


D-TG 


ppt 


170 


10.0 


17.5 


0.7 


4.5 


D-this 


ppt 


170 


7.0 


25.7 


2.9 


4.2 


D'-this 


ppt 


170 


8.6 


22.6 


1.9 


4.9 



Note. — The mass fract ion of heavy elemeii ts is given in solar units Zq = 1.92%. Models 
(A-E)-TG are taken from lGuillot et al.l (120031 ). Present models with prime are designed to 
give the best agreement, see Table [3] for the choice of parameters. 



- 19 - 



Models (A-E)-TG indicate that Ti affects Mc by several Mq, but apparently not sys- 
tematically. Keeping other constraints constant, we always find a slightly higher core mass 
for a larger 1-bar temperature. On the other hand, shifting the layer boundary between the 
envelopes outwards, e.g. from 2.0 to 1.4 Mba r in model A' - this, always magnifies the core 
mass such that the result Mc{A) > Mc{B) of iGuillot et al.l (120031 ) can be reproduced. Our 
value of Mc = 2.8 in model A'-this can be further increased by 0.6 if we use ices as 
core material instead of rocks. 
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Table 3. Parameters of present Jupiter models using SCvH EOS 



Model 


P„, [Mbar] 


Y [%] 


J4/10-^ 


A&B-this 


2.00 


27.5 


-5.84 


E-this 


1.73 


27.5 


-5.84 


D-this 


1.76 


27.5 


-5.84 


A'-this 


1.40 


28.0 


-5.89 


B'-this 


2.00 


27.8 


-5.82 


E'-this 


1.77 


28.5 


-5.79 


D'-this 


1.74 


28.5 


-5.89 



Note. — In models using SCvH-ppt EOS, coincides with the PPT. These models have 
been calculated in order to check the consistency with former results, see Table |2l 
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In agreement with iGuillot et al.l (120031 ) the total abundance of metals Mz = Mz^^^ + 



^^mot + increases by 10-15 when using SCvH-i instead of SCvH-ppt. Since all our 
solutions yield Zmoi < ^met with SCvH-i a r id Zrn n\ > ^met with SCvH-ppt in agreement 
with the bulk of solutions in iGuillot et al.l (120031 ) . the inverse ratio of models A-TG and 
E-TG cannot be reproduced. Apart from the deviations of 43% for Zmoi(E') and 170% for 
Zmct(D')) the deviations of all other solutions are below 30% so that they lie within the 
range of solutions found in IGuillot et~aD J2003h . With < < 10, 17.5 < Mz < 35.3, 
2.3 < Zjnoi < 5.3, and 0.7 < Z^^et < 6.1 the solutions (A-E)-TG span large ranges in units 
of Earth masses and solar abundances Zq = 1.92%, respectively. Our nearest results differ 
by up to 1.4 Mq) for Mc, 5.1 for Mz, 1.0 for Z^oi, and 1.2 Zq for Z^^t from these 
results. To reduce the disagreement in model D with respect to Zmet and Mz, which are 
nearer to model E-TG than to D-TG, the mean helium content should be higher or the EOS 
of metals be more compressible. 

We conclude that our method is able to reproduce former results well enough to smoothly 
continue and complement investigations of the interplay between the internal structure of 
giant planets and the EOS applied. 



4.2. New results with LM-REOS 

Figures [7] and [8] show the full range of solutions using QMD data in comparison with 



models based on SCvH-ppt and SCvH-i adapted from IGuillot et al.l (120031 ) with respect to 
Mc, Mz, Z^oi, and Z^^t- Figures |9] and ITO] illustrate the effect of a varying P^, Ti, and 
J4 on the position of the QMD data based solutions. For all our solutions, Y^^i = 0.238, 
Y = 0.275, Pjn is varied from 3 to 5 Mbar, Ti from 165 to 170 K and J4 from -5.84 to 
-5.89x10"^. Metals are represented either by He4 or by H2O, see §3.51 



The behavior of the isentropes already indicates (see Fig. [6]) that the solutions based on 
LM-REOS are nearer to the SCvH-i than to the SCvH-ppt results. This is a clear consequence 
of the missing PPT in the H-REOS data. The PPT makes the H-EOS less compressible in 
the molecular regime and more compressible in the metallic regime before it is dominated 
by electronic degeneracy. 

The core mass of QMD-based models does not vanish for transition pressures up to 
5 Mbar and J4 < J^^*^. On the other hand, the core mass of our SCvH-i models does 
not exceed 4 M© at transition pressures as low as 1 Mbar. Higher core masses could be 
achieved with SCvH-i for Pm < 1 Mbar, but are not reasonable since the layer boundary 
is to represent the PPT which occurs above 1 Mbar for SCvH-ppt. Higher core masses of 
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QMD-based models can be explained by a smaller compressibility at pressures above 10 Mbar 
as seen from the isentropes. For a given abundance of metals, a smaller compressibility of 
hydrogen adds less mass to the isentrope and leaves more mass for the core. Similary, the 
smaller the compressibility of the EOS of metals (Fig. H]), the higher their abundance needed 
to contribute mass to the isentrope in order to match J2 and J4. 

A pronounced difference between the isentropes calculated with the LM-REOS and 
SCvH-i occurs around 1 Mbar, where the QMD pressure is up to 30% smaller. This feature 
leads to smaller Zmoi and higher transition pressures, see Fig. [HI which can be explained as 
follows. For a given pressure-density relation along the isentrope as required by equations ([9]) 
and ([7j), a higher compressibility of the hydrogen component has to be compensate d by a 



smaller metal abundance in order to keep J2 and J4 at the observed values. Guillot (120051 ) 
has also shown that the mean radius of the contribution of level radii to the gravitational 
moments J2n increases with increasing index n and is placed around the layer boundary 
in case of J2 and J4, whereas the contribution to higher moments than J2 from the deep 
interior is negligible. The molecular region contributes slightly more to J4 than to J2, and the 
metallic layer contributes slightly more to J2 than to J4. If the compressibility of hydrogen 
becomes large near the molecular-metallic layer boundary, Zmoi has to be chosen small in 
order to not exceed | J4I; but simultaneously J2 will become too low and thus Z^^t has to be 
enhanced strongly, again enlarging simultaneously | J4I beyond the desired value. Therefore, 
in order to adjust both gravitational moments with a Z^o\ not too small, the influence of 
Zjnct on J4 has to be reduced by shifting the layer boundary inwards. LM-REOS based 
models of Jupiter have Zmoi = 0.2 — 1.8 Z©, Z^^t = 4.8 — 9.6 Zq for = 3 to 5 Mbar. 



4.3. Systematic behavior of LM-REOS solutions 

In this paragraph we illustrate the influence of the uncertainties in observational and 
free parameters on the position of the solution of LM-REOS based models in the usual M^. vs. 
Mz and Z^oi vs. Zmet diagrams. For our point of reference, displayed as filled circle in Figs. 
|9]and[IOl we have chosen the values 74/10"^ = -5.84, Ti = 170 K, and Pm = 4 Mbar. As 
explained above, Z^oi increases with | J4I. We do not consider solutions for 1741/10*^ < 5.84 
since they result in Zmoi < Zq or Mc = 0. In the molecular envelope where the isentrope is 
more sensitive with respect to temperature than in the degenerate metallic region, a cooler 
interior initiated by a smaller Ti enhances the partial density of the He-H mixture and 
reduces the fraction of metals by 1 M®. Most important, the solutions are affected by the 
choice of the transition pressure. For Pm < 3 Mbar, no solution exists with Zmoi > Zq. 
Solutions with 1 < Zmoi < 1-8 x Zq are found for 3 < Pm < 5 Mbar in combination with 
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Ti = 170 K and J4 < —5.84 x 10~^. A slightly higher enrichment factor of 2.0 is possible 
for Pin = 7 Mbar, before the core mass shrinks to zero. In contrast to the large influence of 
Pjn on Zinoi and ^met, Mz just varies within one M^. Finally, representing metals by H2O 
instead of He4 enhances Z^^t by 50% and shifts Mz from about 30 to about 40 M^. 

4.4. Profiles of chemical components 

We compare the abundances of single chemical species, see §3.71 along the radius inside 
Jupiter as calculated with LM-REOS (model JTTh in Fig. [TT], corresponding to the point of 
reference for H2O in Fig. ^ and with SCvH-ppt (model Jllb in Fig. [TT]). In these figures, a 
full arc segment corresponds to 100% in mass. The radius coordinate is not displayed here 
but scales linearly from the center to the pressure level of 1 bar. In both of the models the 
core consists of rocks, 74/10""^ = —5.84, Y = 27.5%, and Ti = 170 K. There is no degree of 
freedom in model Jllb and the only degree of freedom in model Jllb. P^, is set to 4 Mbar 
in order to give Z^^i > 1- The three most important differences are: (i) the steep onset of 
dissociation and ionization of H-atoms with the PPT in model JTTb wheras LM-REOS yields 
a smooth increase of dissociation; (ii) a high transition pressure of 4 Mbar in model JTTk. 
well beyond the location of the PPT in SCvH-ppt; (iii) a ratio of heavy elements Z^oi'-Zmet 
of 4:3 in model Jllb compared to 1:6 with LM-REOS. It is interesting to note that once 
neutral H and He"*" are formed in the SCvH-ppt model, they are almost immediately ionized 
further, whereas the fraction of ionized He of only 6% in the deep interior in model JTTh is 
not resolved. 

For model Jllb. an ASCII data table containing the profiles of pressure, temperature, 
density, composition and the figure functions along the radius is availiable in the electronic 
edition of this Journal. A shortened version of this data table is shown in Tab. HI where the 
five rows, from top to bottom, present the 1-bar-level surface, the transition from the outer 
to the inner envelope, and the layer transition from the inner envelope to the core. 



Table 4. Jupiter model Jlla 



m 
[M®] 


P 

[Mbar] 


I 


T 
[K] 


P 

[gcm~^] 


Y 

[%] 


Z 

[%] 




Xhc 




S2 


Si 


S6 


317.8336 


l.OOOOE-06 


1.00000 


1.7000E+02 


1.6653E-04 


23.307 


2.072 


1.0000 


1.0000 


0.0000 


-4.501E-02 


1.984E-03 


-2.470E-04 


223.6070 
223.5810 


3.9998E+00 
4.0013E+00 


0.72384 
0.72380 


8.8682E+03 
8.8694E+03 


1.3239E+00 
1.5253E+00 


23.307 
24.476 


2.072 
16.616 


0.0095 
0.0095 


1.0000 
1.0000 


0.0000 
0.0000 


-3.522E-02 
-3.522E-02 


1.062E-03 
1.062E-03 


-1.742E-04 
-1.744E-04 


2.75441 
2.75308 


3.8385E+01 
3.8393E+01 


0.08426 
0.08426 


1.8571E+04 
1.8572E+04 


4.3293E+00 
1.8037E+01 


24.476 
0.000 


16.616 
0.000 


0.0000 
0.0000 


0.9405 
0.0000 


0.0574 
0.0000 


-1.160E-02 
-1.160E-02 


1.063E-04 
1.063E-04 


-7.278E-06 
-7.278E-06 



Note. — This table is a truncated version of a machine-readable table that is published in its entirety in the electronic edition of the Astrophysical Journal. 
Column headings from left to right: mass coordinate, pressure, level coordinate (radius, see Eq. [S]!, temperature, mass density, mass fraction of helium, mass 
fraction of metals (H2O), particle fraction of H2 molecules with respect to the H subsystem, particle fraction of neutral He (column 9) and of singly ionized 
He (column 10) with respect to the He subsystem, and columns 11-13: the dimensionless figure functions S2, S4, sg. 
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5. Evaluation of the new results 

In this chapter we discuss to what extent the new results for Mc, ^met, ^moi, and Pm 
obtained by using LM-REOS are in agreement with experimental EOS data (e.g. principal 
Hugoniot curve), evolution theory, element abundances, and H/He phase separation. 



5.1. Core mass Mr 



Saumon & Guillot (120041 ) found that H-EOS with a small maximum compression ratio 
of only 4 along the Hugoniot curve can yield small core masses lower than 3 M® and com- 
ment that the apparent relation between the stiffness along the principal Hugoniot and the 
core mass may be not unique. In the same sense, our H-REOS reproduces the maximum 
compression ratio of 4.25 as derived from shock-wave experiments, but the core masses range 
up to 7 M®. 

To understand the indirect effect of the compressibility k, of hydrogen on the core mass, 
we study in particular its impact in the deep interior (/tmet), around the layer boundary 
between the envelopes (Km), and in the outer molecular region (/tmoi)- The core mass depends 
essentially and directly on the mass density pmet in the deep interior: The higher Pmet, the 
lower Mc. Clearly, pmet can either be enhanced by Kmet or by Z^et- For example, the case 
of a smaller Kmet leading to a smaller pmet occurs with LM-REOS. Furthermore, Z^et is 
diminished by both k,^, which reduces the need for metals in order to adjust J2, and the 
Zmoi chosen to reproduce J4. Finally, Z^oi decreases with Kmoi, for instance in case of SCvH- 
ppt. Due to this propagation of effects, the behavior found in ISaumon fc GuillotI (120041 ) 
may cor respond to the coinciden ce of a small Km an d a large /tmet. betw e en 5- 15 Mbar, see 
Fig. 1 in I Saumon fc GuillotI (l2004j ). In agreement with lSaumon fc GuillotI (120041 ) we conclude 
that the compressibility along the principal Hugoniot curve, which is restricted to densities 
below about 1 g cm~^, does not determine the core mass alone. Experimental data for the 
hydrogen EOS off the principal Hugoniot curve, i.e. near the isentrope, are in this context 
urgently needed. For this, new experimental techniques such as reveberating shock waves or 
precompressed targets can be applied. 



5.2. Abundance of metals Z^et 



In our new Jupiter models with LM-REOS, Zmet is enriched over solar abundance by 
a factor of 5 to 10 and exceeds Zmoi by a factor of 4 to 30. This featur e is consistent with 
the standard giant planet formation scenario, the core accretion model (lAlibert et al.ll2005l : 
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Pollack et al.lll996l ). where the planet grows first by accretion of planetesimals onto a solid 



core embryo. If the core has grown such that surrounding nebula gas is attracted, an envelope 
forms an d the planet grows b y accretion of both gas and planetesimals, either sinking towards 
the core (jPoUack et al. 1996 ) or dissolving in the envelope. If the luminosity of the envelopes 
has reached a critical value, the energy loss due to radiation cannot be supplied anymore by 
infalling planetesimals, and the whole planet starts to contract with an even enhanced gas 
ac cretion rate (run-aw ay growth). Depending on a variety of parameters described in detail 



m 



Alibert et al.l (120051 ). at the end of the lifetime of a protoplanetary disk a Jupiter-mass 



giant planet may have formed with a core mass of several and a total mass of heavy 
elements of 30-50 My. Further evolution of t he planet may include erosion of core material 
( jChabrier &: Baraffd 120071 : iGuillot et al.ll2003[) an d homogeneous redistribution of destroyed 
planetesimals due to convection (iKerleyl l2004al ). Both of these processes depend on the 
ability of convection to overcome compositiorial gradients and the ir efficiency is not known 
in detail up to now (jChabrier fc Baraffd 120071 : iGuillot et al.l 120031 ). Thus, in the framework 
of the core accretion model, the distribution of heavy elements in our Jupiter models can 
qualitatively be explained as a consequence of planetesimal deliveries in the deep interior 
and their redistribution within the metallic layer by convection. 



5.3. Abundance of metals Z 



mol 



With the exception of oxygen, measurements of the abundance of the elements C, N, S, P 
and the nob l e gase s beside He and Ne give two- to fourfold enrichment over solar abundance 



see 



GuillotI (120051 ) for an overview). On the other hand, the Galileo probe gave an 0/H 



ratio of only 0.2 x Zq. It has been argued that the probe fell into a non-representative dry 
region and the O abundance was still rising with depth when the probe stopped working 
( lKerleyll2004al ). Our result of Zmoi = 1 — 1.8 x Zq is too low by a factor of two. We note 
that LM-REOS based three-layer Jupiter models are difficult to reconcile with the observed 
atmospheric abundances of heavy elements. 



5.4. Layer boundary Pm 

None of the mechanisms known to affect e.g. H/He phase separation or the prop- 
erties of the nebula in the neighbourhood of the forming planet, could explain a depletion of 
heavy elements from the gas phase below solar abundance. To let the Jupiter models have 
Zmoi > Zq requires Pm > 3 Mbar using LM-REOS. In the three-layer model the discontinuity 
of the fraction of metals is assumed to occur at the same pressure as the discontinuity of the 
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fraction of helium. For helium, the discontinuity has be en argued to coincide with a PPT 
of hydrogen or with the onset of H/He phase separation (IStevenson &: Salpeterlll977bl ). Our 
QMD simulations show no indications of a PPT so far. 

H/He phase separation on the other hand is a consequence of non-linear mixing and 
characterized by a line T{x, P) of demixing which denotes for a certain concentration x 
the maximal, pressure-dependent temperature for which the system shows phase separation. 
While neutral He atoms are expected to separate from metallic hydrogen, ionized helium 
becomes miscible again at high pressures. However, only few theoretical re sults for the 
demixing line exist up to now. Latest simulations (jPfaffenzeller et al.l Il995l ) have a line 
of demixing at temperatures well below the adiabatic temperature profiles of Jupiter and 
even of Saturn. Some theoretical EOS of H/He mixtures predict that t he Jupiter isentrope 
intersects the high pressure boundary of the dern ixing line between 1.5 (IFortney fc Hubbard 
20041 ) and 7 Mbar (IStevenson fc Salpeteiill977al ). At smaller and at higher pressures along 
the isentrope, the system would mix again. These results are not confident enough to support 
or exclude a layer boundary around 4 Mbar. 



5.5. Possible improvements 

Three-layer models of Jupiter based on QMD data are consistent with the core accretion 
model, but cannot reproduce observed heavy element enrichments by a factor of 2 to 4. 
Transition pressures below 3 Mbar are even accompanied by Z^oi < I Zq. Considering the 
agreement of our ab initio QMD results with measured Hugoniot curves and reflectivities 
for hydrogen, the present results for Jupiter point to the importance of three topics to be 
addressed in the future: (i) inclusion of mixing effects in the EOS instead of applying linear 
mixing, (ii) recalculation of the H/He demixing line within an ab initio approach; and (iii) 
adjustment of the structure model in order to account for phase boundaries of H/He demixing 
and He enrichment at deeper levels while still keeping the model as simple as possible. 

First r esults on the eff e ct of nonlinear mixing to the pressure-density relation were 



obtained by lVorberger et al.l (120061 ) within DFT-MD simulations. They find that the volume 
at constant pressure is enhanced up to 5% at temperatures typical for the deep molecular 
layer where our H-EOS exhibits the highest compressibility, i.e. the smallest volume relative 
to SCvH-i. Taking into account the small mass fraction of metals, a compensation of a 2% 
density reduction in average of the H/He component in the outer envelope would require 
the fraction of metals to be doubled to about 3 Zq which is just the average of observed 
particle species. This effect is, therefore, a serious candidate to remove the mismatch between 
calculated and observed Z^^i, see discussion above. 
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Using LM-REOS, the NMT in Jupiter occurs at 0.5 Mbar and is clearly separated from 
the layer boundary around 4 Mbar with about 80 in between. If the order of the latter 
one will be confirmed by future work on the H/He phase diagram, the possibility of an 
extended inhomogenous layer enhancing the compressibility should be taken into account in 
the structure model. 



In order to complement these and former investigations (e.g. Hubbardlll969 : Saumon fc Guillot 



2004 ) of the interplay between the internal structure of giant planets and the EOS applied, 



a calculation of the evolution history including H/He phase separation is also necessary. 



5.6. Impact for other planets 



Corresponding results can be obtained also for Saturn, Uranus, and Neptune but with 
higher uncertainty regarding Mc, Zj^^i and Zj^^t- For most of the extrasolar giant planets 
(EGP) detected so far the only structural parameter constrained at all is the minimum 
mass. Transit detectio ns also allow to determine the mean radius with an error of about 
10% (ISato et al.l |2005| ) resulting mainly from the uncertainty in the stellar radius. Even 



the few transiting planets detected so far have revealed a huge variety of mass-radius rela- 
tionships. To explain them and to derive implications for their formation and contraction 
history, well-grounded knowledge of the Solar system giant planets is crucial. For instance, 
characteri stics as helium depletion and compositional gradients are expected t o apply also 
to EGPs (IChabrier fc Baraffel[2007l : IPortnev fc Hubbard! [20oi : ISato et all 120051 ). 



6. Conclusions 

We have computed three-layer structure models of Jupiter using new EOS data for H, 
He, and H2O. They are composed of ab initio data in the warm dense matter region covering 
at least 97% of Jupiter's total mass and of chemical model EOS at low densities. With 
to 7 M0, the core mass compares with previous results satisfying the same observational 
constraints. Most of the other properties of QMD-based models such as (i) up to tenfold 
enrichment with heavy elements in the metallic layer, (ii) small heavy element abundance in 
the outer envelope, and (iii) a transition pressure of at least 3 Mbar, are a consequence of 
the high compressibility of the hydrogen EOS between 0.1 and 10 Mbar along the isentrope, 
but not along the Hugoniot. Transition pressures of 3 to 5 Mbar agree with estimates of the 
onset of remixing of hydrogen and helium that are demixing at smaller pressures. 

These results clearly underline the importance of calculating the H/He phase diagram 
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with respect to the EOS and the region of hehum immiscibihty in order to improve our 
understanding of the internal structure of hydrogen-rich giant planets. These issues will be 
addressed in future work. 

We thank J. Kohn, W. Lorenzen, T.R. Mattsson, I. losilevskiy, and T. Guillot for 
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Ph [g cm"^] 

Fig. 1. — (Color online) Isotherms for hydrogen for different EOS tables: FVT"*", QMD, 
and H-REOS. The thick line shows the pressure-density relation of the hydrogen fraction 
of a typical Jupiter isentrope and the grey area masks the non-accessible region below the 
cold-curve represented here by a 100 K isotherm of QMD data. 
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Fig. 2. — (Color online) Isotherms for helium for T=4000, 15800, 31600 K and three different 
EOS tables: Sesame 5761, QMD data, and He-REOS. The thick line shows the pressure- 
density relation of the helium fraction along a typical Jupiter isentrope. 
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Fig. 3.— (Color online) Isotherms for H2O for T=200, 500, 1000, 4000, 6000, and 20000 K 
according to different EOS tables: a combination of QMD and Sesame 7150 for T > 100 K, 
a combination of WP and Sesame 7150 for T < 1000 K, FW iFeistel fc Wagneil kood ) for 
ice I, Sesame 7150, and QMD data. The thick line shows the pressure-density relation of 
the H2O fraction along a typical Jupiter isentrope. All these isotherms are incorporated into 
H2O-REOS. 
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Fig. 4. — (Color online) Isotherms of the Z-component for T=4000, 6000, and 20000 K 
according to different EOS: Sesame 7150 for H2O, QMD data for H2O, H2O-REOS, and 
He4-RE0S. 
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Fig. 5. — (Color online) Coupling parameter F and degeneracy parameter G inside Jupiter 
for different EOS: LM-REOS, SCvH-i, and SCvH-ppt. 
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Fig. 6. — (Color online) Hydrogen adiabats for Jupiter, determined by T=165 K at P=l bar, 
computed with three different hydrogen EOS: H-REOS, SCvH-i, and SCvH-ppt. The scale 
on the right shows absolute pressures and the scale on the left shows relative differences in 
pressure with respect to the SCvH-i-adiabat. 
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Fig. 7. — Core mass and total mass Mz of heavy elements in units of Earth masses for 
three different EOS. Solutions with LM-REOS (solid lines) demonstrate a strong influence 
of the choice of the EOS of metals to the resulting total heavy element abundance and are 
displayed in seperate boxes labeled 'He4' and "H2O', respectively. 
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Fig. 8. — Abundance of heavy elements in solar units in the molecular layer (Zmoi) and in 
the metallic layer (Zmet) for three different EOS. The dotted line indicates equal abundances. 
Solutions with LM-REOS (solid lines) demonstrate a strong influence of the choice of the 
EOS of metals to the resulting value of Z^et and are displayed in seperate boxes labeled 
'He4' and 'H2O', respectively. 
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Fig. 9. — Same as Fig. [7] but only for solutions using LM-REOS. The arrows indicate the 
shifts of the reference solution (Ti = 170 K, Y = 0.275, Pm = 4 Mbar, J4/10^ = -5.84) if 
Ti is decreased to 165 K, | J4I increased by 1 a, or enhanced from 3 to 5 Mbar. 
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Fig. 10. — Same as Fig.[H]but only for solutions using LM-REOS. See Fig.[9]for a description. 
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Fig. 11. — Schemes of Jupiter models satisfying the same constraints; left: model JTTk. 
right: model JTTb. At the layer boundaries the values of pressure and temperature are 
given. The abundances of metals and of chemical species along the radius are indicated by 
grey scales. An arc segment corresponds to 100% in mass. For model Jllb. an ASCII data 
table containing the profiles of pressure, temperature, density, composition and the figure 
functions along the radius can be found in the electronic edition of this Journal. 



